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SUMMARY 

We introduce the Bozzolo-Ferrante-Smith (BFS) method for alloys as a computationally effi- 
cient tool for aiding in the process of alloy design. An intuitive description of the BFS method is 
provided, followed by a formal discussion of its implementation. The method is applied to the 
study of the defect structure of NiAl binary alloys. The groundwork is laid for a detailed progres- 
sion to higher order NiAl-based alloys linking theoretical calculations and computer simulations 
based on the BFS method and experimental work validating each step of the alloy design process. 
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INTRODUCTION 


A number of industries including the electronics, computer, aeronautics, aerospace, and other 
transportation related enterprises have long relied on the application of new high-performance 
materials for their competitive advantages. Yet, in aeronautics and aerospace applications in par- 
ticular, the potential loss in lives and liability that could result from any material failure dictates 
that these systems must be very well characterized and understood. Therefore, not only do these 
new materials have to push the envelope in terms of performance, at the same time they need a 
high degree of maturity. This makes the introduction of new materials in these industries time 
consuming and costly due to the need for painstaking development, testing, and characterization. 
At the same time, without the large defense-related expenditures that were possible during the 
cold war era, it is necessary to find new ways for industry to operate that will provide significant 
reductions in expense and cycle time in the development, maturation, and certification processes 
for new materials. To achieve these goals, the experimental approach to materials research and 
design will have to be augmented with computational processes. The goal is to develop computa- 
tional capabilities for alloy design with crystalline materials to the same or greater degree than 
exists today in the polymer and pharmaceutical industries. 

In a series of articles we will describe the integration of a new computational materials model, 
known as the BFS method [1], into an actual materials design program. The intent of this inte- 
grated experimental and computational program is to develop revolutionary new alloys based on 
the ordered intermetallic compound NiAl as a replacement for Ni-based superalloys in aero-tur- 
bine applications. NiAl offers a number of advantages over Ni-based superalloys, including a 
higher melting temperature, better environmental resistance, significantly higher thermal conduc- 
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tivity, a better response to thermal barrier coatings, and a reduction in density by as much as a 
third without significantly impacting creep strength [2], However, the development of an alloy 
composition with the correct balance of high temperature strength and low temperature fracture 
resistance has continued to elude investigators [3], leading to continued alloy design-structure- 
property investigations and the eventual integration of computational processes (BFS method) in 
the search for an acceptable alloy. 

Through the eventual partnership between experimental and theoretical procedures, we were 
able to determine the naturally evolving role of computational methods in an alloy design pro- 
gram. This role is schematically illustrated in Figure 1 . At first, the design approach reported in 
this work was guided by past experimental work of an empirical nature. For most alloy design 
programs this is essentially one continuous loop (Figure l.a). But an integrated alloy design pro- 
gram that augments experimental results with theoretical modeling would follow the path 
described in Figure l.b. We started by modeling what was known about the NiAl system. This 
allowed us to establish confidence in the model, make changes and optimize the model, and to 
help understand and interpret the results. The model also provided the energetic explanations for 
many of the experimental observations, thus giving us greater insight into the alloy system. As 
confidence in the model grew, we started filling the gaps in our experimental database with com- 
putational results. Finally, as alloy compositions became more complicated and burdensome to 
evaluate, we used the computational model to direct the experimental work, and only experimen- 
tally verified those few compositions that looked promising or that we needed in order to refine 
the model. The advantage of this integrated approach is that it funnels most of the alloys through a 
virtual production and screening process significantly reducing the number of alloys actually pro- 
duced and characterized. This saves significant time and expense since the experimental process is 
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Figure l.b Integrated Alloy Design Program with Mutual Feedback 
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the rate and resource limiting step. 

This paper represents the first in the series of articles detailing our design program on NiAl- 
based alloys. In forthcoming papers, we will describe our experimental and theoretical investiga- 
tion of a series of 3-, 4- and 5-element alloys based on NiAl, highlighting the connection between 
the two approaches. In this first paper we will provide a detailed description of the BFS method 
and use it to reproduce the basic features of the base material, with as much accuracy as possible, 
since the forthcoming extension of the model to complex systems will rely only on our ability to 
describe the most basic system. 

Therefore we will initially concentrate on the P-phase of the binary Ni-Al system, which 
exists over a range of stoichiometry from about 45-60 at. % Ni [4], At the stoichiometric compo- 
sition NiAl should exist in a perfectly ordered state where the Ni atoms occupy the cube comers 
and the A1 atoms occupy the cube centers of a generalized body-centered cubic lattice. Ni-rich 
alloys are characterized by antisite point defects, where Ni atoms occupy sites in the A1 sublattice, 
resulting in a decrease in lattice parameter and an increase in density with increasing Ni content. 
A different behavior is observed on the Al-rich side of stoichiometry. There is a steep decrease in 
lattice parameter as well as a greater than expected decrease in density with increasing A1 content. 
The presence of vacancies in Ni sites would explain such behavior [4-8]. Recent X-ray diffraction 
experiments [9] suggest a richer structure: the evidence, while strongly favoring the presence of 
vacancies in Ni sites, also suggests the possibility of some vacancies in A1 sites in a 3:1 ratio. 
Moreover, local ordering of vacancies may be preferred over a random distribution of individual 
point defects [6]. 

While some of the concepts introduced in this first paper might seem too detailed for the task 
at hand (i.e. describing the defect structure of P-NiAl), they will be strictly necessary for a full 
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understanding of the results from more complex systems that will appear in subsequent papers in 
this series. In other words, we are establishing the validity of our model by eliciting a positive 
response to ‘Q2’ in Fig. l.b before proceeding further along the flow chart. 

In most cases, theoretical techniques at the atomic level are used for validating or understand- 
ing experimental findings, a task somewhat removed from the actual process of alloy design. In 
this series, we embark on a process of materials development starting with the validation of our 
theoretical approach with established experimental knowledge of the basic, binary system (this 
first paper), building on our theoretical predictions as the complexity of the system increases (ter- 
nary, quaternary alloys), finally allowing us to determine the structural properties of a 5-element 
alloy (last paper) based on this groundwork. As a whole, this series of papers is meant to prove 
another point, which goes beyond the features of the particular problems at hand: that semiempir- 
ical methods, tied to powerful numerical techniques and implemented in the framework of an eco- 
nomical and easy to handle operational procedure, can constitute a very useful tool in effectively 
aiding in the area of materials design. 

This paper introduces, in detail, the BFS method and its operational equations. After introduc- 
ing the novel way of partitioning the alloy formation process proposed by the BFS method and 
conceptually defining the quantities of relevance, we present details of the operational procedure. 
We continue with an example illustrating the implementation of the method. As a consequence of 
the tremendous advance in semiempirical techniques over the last decade, we conclude this sec- 
tion with a brief description of recent work done with other techniques of comparable complexity 
as well as their relation to BFS. We continue by introducing the computational scheme chosen to 
address the structure of NiAl, a procedure which will be used and expanded in subsequent papers 
dealing with higher order alloys. We conclude by providing the BFS predictions for the defect 
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structure and physical properties of (3-NiAl, comparing them with experimental results when pos- 
sible, thus providing us with confidence in our model. 


THE BFS METHOD 


The BFS method is one of the latest additions to the growing family of quantum-mechanical, 
approximate techniques for the calculation of materials properties. It has a strong foundation in 
quantum theory, but simplifies the numerical effort usually associated with ab-initio methods. 
While a great deal of detail (i.e. electronic structure) is lost when applying semiempirical methods 
to specific systems, valuable insight on the general behavior of physical systems is gained, since 
the simplicity of the formulation allows for a quick and sufficiently accurate estimate of general 
properties. 

Several methods successfully deal with single -element metallic systems, providing a wealth of 
information on their physical properties and defect structures [10-15J. However, in the case of 
extended defects, i.e. surfaces, some methods exhibit limitations in their predictive power [11]. In 
most cases the corresponding potentials or the parameters used in the algorithm, are usually deter- 
mined from bulk properties and are therefore limited in their ability to properly reproduce these 
defects. However, the overall power of semiempirical techniques has enabled a better understand- 
ing of defects such as surface phenomena in spite of these shortcomings, which has lead to exten- 
sive growth in the area of computer simulations of defects. 

For alloy structures, semiempirical methods have been less effective in providing results with 
the same degree of reliability as those for single-element systems. This is partly due to the nature 
of the formalisms. Either an impractical number of parameters is needed (modified embedded 
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atom method [11] and tight-binding methods [12], for example) or specific potentials have to be 
determined on a case-by-case basis to deal with each specific alloy composition or phase, there- 
fore limiting their transferrability. 

In spite of these limitations, there has been remarkable progress in the last few years, hi com- 
bination with the development of powerful first-principles methods and an ever increasing compu- 
tational power, theoretical approaches are becoming tools of great value for understanding, 
modeling and designing materials suited for industrial applications. 

Perhaps the main contribution of the BFS method to the field of computationally oriented 
semiempirical methods, is that it is basically free of most of the restrictions that limit the applica- 
tion of comparable techniques [16]: there is no restriction on crystal structure or the type of 
atomic species considered and the number of input parameters is reduced to a minimum. More- 
over, no experimental input is necessary, as all the input parameters can be determined by simple 
first-principles calculations [17]. 

These advantages of the BFS method make it particularly suitable for broad alloy design prob- 
lems where the nature of the predictions is comparable to the information that can be obtained 
from experimental analysis (in this work, for example, we attempt to study the structure of 5-ele- 
ment alloys, both theoretically and experimentally). It is precisely this issue that motivates the 
present work: to develop a reasonably accurate and computationally simple theoretical approach 
that can provide valuable input for the development of new materials, in close conjunction with 
experiment. 

Since its inception a few years ago [18], the BFS method has been successfully applied to a 
variety of situations ranging from segregation profiles [16], alloy structure [19], surface alloying 
of immiscible metals [20], to numerical simulations of the scanning tunneling microscopy tip- 
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sample interaction [21], bulk alloy design [1,22] and alloy surface structure [23]. Two of the char- 
acteristics of the early applications of the BFS method, i.e., the need for experimentally deter- 
mined input data and complete reliance on the pure constituent properties [18], initially imposed a 
limitation when dealing with systems like P-NiAl, where both constituents are fee elements but 
the alloy is bcc-based. No appropriate experimental input for bcc-Ni and A1 is available. There- 
fore, as shown in detail later in this work, we have reformulated the method relying only on pure 
first-principles determined input, thus freeing BFS from limitations imposed by the restricted 
availability of experimental data. We also avoid the potential problem of inconsistency or inaccu- 
rate data obtamed from different experimental techniques [14]. While the method retains the same 
operational algorithm used hi its previous applications, the new parameterization scheme pre- 
sented in this paper greatly enhances its range of applications and accuracy. 

General concepts 

The BFS method provides a simple algorithm for the calculation of the energy of formation 
A H of an arbitrary alloy (i.e, the difference between the energy of the alloy and that of its individ- 
ual constituents). In BFS, the energy of formation is written as the superposition of individual 
contributions of all the atoms in the alloy. 


AH = S< £ ':- £ :> = 2A (1) 

I I 

where E j ’ is the energy of atom i in the alloy and E i its corresponding value in a pure equilibrium 
monatomic crystal. In principle, the calculation of A H would simply imply computing the energy 
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of each atom in its equilibrium pure crystal and then its energy in the alloy. In BFS, beyond 
directly computing the difference e, for each atom in the alloy, we introduce a two-step approach 
for such a calculation in order to identify contributions to the energy due to structural and compo- 
sitional effects. Therefore, we break up the individual contributions of each atom e f to the total 
energy of formation AH of the sample into two components: a strain energy and a chemical 
energy. 

While there is a certain level of arbitrariness in how this separation is implemented, it is only 
meaningful when a good representation of the actual process is obtained by properly linking both 
contributions, hi other words, a proper coupling of these two apparently independent processes 
(strain and chemical effects) must be accomplished in order for the final result to approach the 
result one would obtain if a straightforward calculation (i.e. ab-initio methods) was carried out. 

We define the BFS strain energy as the contribution to the energy of formation from an atom 
in an alloy computed as if all the surrounding atoms were of the same atomic species, while main- 


(a) 


(b) 


(c) 




Fig 2: (a) A pure, equilibrium crystal (reference atom denoted by the arrow), (b) a reference atom (denoted 
by the arrow) in the alloy to be studied (atoms of other species denoted with other shading) and (c) the 
same reference atom in a monatomic crystal, with the identical structure of the alloy to be studied, but with 
all the atoms of the same atomic species as the reference atom, for the calculation of the strain energy tenn 
for the reference atom. The strain energy is the difference in energy of the reference atom between (c) and 
(a). 
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taining the original structure of the alloy. To visualize this concept, Fig. 2. a represents the atom in 
question (identified with an arrow) in an equilibrium position in its ground state crystal (arbi- 
trarily represented by a simple cubic lattice). Fig. 2.b shows the same atom in the alloy being stud- 
ied (also arbitrarily represented by a different crystallographic symmetry). Two things can be 
different between the reference crystal and our alloy. First, atoms of other species may occupy 
neighboring sites in the crystal and second, the crystal lattice may not be equivalent in size or 
structure to that of the ground state crystal of the reference atom. In Fig. 2.b, the different atomic 
species are denoted with different symbols from that used for the reference atom, and the differ- 
ences in size and/or structure are denoted with a schematically different atomic distribution as 
compared to the ground state crystal shown in Fig. 2. a. The BFS strain energy term accounts for 
the change in energy due only to the change in the geometrical environment of the crystal lattice 
(from 2. a to 2.b), ignoring the additional degree of freedom introduced by the varying atomic spe- 
cies in the alloy. In this context, Figure 2.c shows the environment ‘seen’ by the reference atom 
when computing its BFS strain energy contribution. The neighboring atoms conserve the sites in 
the actual alloy (fig. 2.b), but their chemical identity has changed to that of the reference atom 
(fig. 2. a), thus simplifying the calculation to that of a single-element crystal. The BFS strain 
energy term represents the change in energy of the reference atom in going from the configuration 
denoted in Fig. 2. a to Fig. 2.c. 

This choice for the BFS strain energy introduces two main advantages. First, by transforming 
(for each atom) the actual alloy into a monatomic crystal, it greatly simplifies the calculation of 
the energy of the reference atom in the alloy structure, making it amenable to a large number of 
theoretical techniques that can efficiently deal with this simplified situation. Second, it gives par- 
tial information concerning only the structure of the alloy, which could serve later to identify fine 
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geometrical effects on structure. 

The second contribution to the BFS energy of formation is the chemical energy. Here, we are 
interested in isolating the effect of compositionally different atoms occupying neighboring sites to 
the reference atom. In order to do this, we leave out any structural information from the original 
lattice by forcing the neighboring atoms to occupy equilibrium lattice sites corresponding to sites 
in the pure cell of the reference atom, changing the composition of the atoms to match the chem- 
ical profile in the alloy lattice. Fig. 3. a shows the reference atom in the actual alloy (similar to Fig. 
2.b), while Fig. 3.b indicates the atomic distribution used in computing the BFS chemical energy 
(note that the lattice used in Fig. 3.b corresponds to that of the ground state crystal of the reference 
atom as shown in Fig. 2. a). 

The chemical energy is then based on the difference between the energy of the reference atom 
in Fig. 3.b and its energy in its ground state crystal (Fig. 2. a). In order to completely free the 
chemical energy from structural defects, certain provisions have to be made which will be clearly 


(a) 



(b) 



Fig. 3: (a) The reference atom (denoted 
by an arrow) in the actual alloy environ- 
ment and (b) the reference atom sur- 
rounded by a chemical environment 
equivalent to that in (a) but with the dif- 
ferent neighboring atoms occupying 
equilibrium lattice sites corresponding to 
the ground state of the reference atom. 
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detailed later in this section. 


So far, we have exploited the flexibility in partitioning the energy of formation to provide 
ourselves with a simple and straightforward framework for computing both the BFS strain and the 
BFS chemical energy contributions. However, a simple addition of these two quantities does not 
necessarily constitute a good approximation to the actual energy of formation. After decoupling 
the structural and chemical processes involved in alloy formation, we now recouple them by intro- 
ducing a coupling function which, in its definition, has the necessary ingredients to ensure that 
this partitioning scheme accurately reproduces a full calculation of the energy of formation when 
recombined in this manner. The definition of the coupling function will be clear once we intro- 
duce the necessary tools for computing the individual (strain and chemical) contributions in the 
following subsections. We therefore write the energy of formation of the alloy as 


Aff = = 5>s + * e c> (2) 

where e s is the strain energy , e c is the chemical energy, g is the coupling function and the sum 
extends over every atom in the alloy. 

Strain energy 

The actual calculation of the BFS strain energy for each atom is straightforward and any tech- 
nique designed to compute the energetics of pure crystals should be appropriate. However, in spite 
of the apparent simplicity of this task, few techniques are actually capable of performing it suc- 
cessfully in general situations. The choice should be limited to only those techniques that intro- 
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duce a substantial level of accuracy, regardless of the geometry of the environment, as well as 
those characterized by straighforward numerical determination. While first-principles calculations 
would be ideal for performing this task, they are still limited by the excessive computational effort 
associated with complex geometries. It is true that because of the simplicity of the NiAl case, 
first-principles techniques could be used without major difficulties. However, this might not be 
true when multiple alloying additions are considered. On the other hand, the computational 
effort involved is not an issue with other semiempirical techniques comparable to BFS, but most 
still fail in properly dealing with this situation in the presence of extended defects from a physical 
standpoint or by relying heavily on parameters that might not be applicable or effective for the 
type of alloy structure under consideration. 

hi all previous applications of BFS, we have used the Equivalent Crystal Theory (ECT) 
[24,25], which in a straightforward manner, provides a simple algorithm for the calculation of 
defect energies even in the presence of such radical defects as surfaces. The underlying concept 
of ECT is simple and clear and it becomes particularly easy to understand if presented in terms of 
the Universal Binding Energy Relation (UBER) [26]. Consider the ground state crystal, character- 
ized by the equilibrium value of the Wigner-Seitz cell, r WSE . An isotropic expansion or compres- 
sion of such a crystal from its equilibrium shape (where all bonds are expanded or compressed 
equally in all directions) yields a universal (in the sense that it has the same shape for all ele- 
ments) curve for the energy vs. Wigner-Seitz radius (or atomic volume), as shown in Fig. 4. a. 

If the reference atom is in an equilibrium position in the ground state crystal, it sits at the bot- 
tom of this curve, with an energy -E c (cohesive energy). If this crystal is deformed anisotropically, 
or if a defect is introduced hi the vicinity of the reference atom (i.e., a vacancy, an interstitial 
atom, etc.), it raises its energy. The curve indicates that there are two values of the Wigner-Seitz 
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radius for which each atom in a homogenous, ideal, crystal have the same energy as the reference 
atom in the defect crystal. We will call these two crystals, equivalent crystals of the reference 
atom. These two positions represent two different conditions. To the right of the minimum the 
crystal has an increased lattice parameter and thus a reduced electron density. To the left of the 
minimum the opposite is true. 


(a) 

Ground state 





Fig. 4: Universal binding energy relation (E vs. r) for a monatomic crystal. When all atoms are in their 
equilibrium lattice sites, the energy of each atom is the cohesive energy, indicated with the arrow in (a). In 
the presence of a surface (b) or in the presence of an interstitial atom (c), a given atom will have a higher 
energy along the binding energy curve, which also corresponds to a certain ideal perfect crystal with an 
expanded (b) or compressed (c) lattice parameter. Any atom in these equivalent crystals have the same 
energy as those of the reference atoms (denoted by an arrow) in (b) and (c). 
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For example, an atom on the surface of a metal (Fig. 4.b), after losing a large number of neigh- 
bors, thus increasing its energy substantially due to the decrease in electron density, has the same 
energy in that situation that it would have hi a perfect crystal but with a larger lattice parameter. In 
terms of its energy, the atom would not distinguish between being on the surface or being in a 
larger version of the ground state crystal. Conversely, if an interstitial atom is introduced close to 
the reference atom, its energy increases to an amount equivalent to the one that it would have in a 
properly homogenously compressed ideal crystal. Figure 5 reiterates these concepts, showing the 
equivalent crystals of a surface atom. As a consequence, every point along the UBER of a certain 
crystal is degenerate, as there is a large number of defects that would have the same energy, and 
therefore, the same equivalent crystal. In other words, the reference atom can find itself in the 
presence of a number of different defects that raise its energy by the same amount, therefore 
assigning the same equivalent crystal to all those defects. 



plane 


|^-ct()+X-*>-| 



Atom in the bulk of an expanded 
(equivalent) crystal 


Figure 5: Equivalent crystals of a surface atom. The surface atom is in an environment of reduced electron 
density, equivalent to that found in an isotropically expanded bulk crystal. 
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ECT provides an efficient algorithm for finding the equivalent crystal for each atom, corre- 
sponding to a certain environment, based on quantum perturbation theory. We refer the reader to 
Ref. 24 for a derivation of the ECT equations from perturbation theory, and to Ref. 25 for more 
details on the implementation of the ECT operational equations. Once the reference atom is cho- 
sen, all that is needed is to know the distances between the reference atoms and a group of its 
neighbors. Although no preset cut-off distance is defined in ECT, it is sufficient in most cases to 
consider all neighbors up to second or sometimes third-neighbor distances. Finding the equiva- 
lent crystal for a certain atom near a defect translates into finding the equivalent lattice parameter, 
i.e., the expanded or compressed lattice parameter of the equivalent crystal, such that the energy 
of an atom in that crystal is the same as that of the atom in the defect crystal. This is done by solv- 
ing the following transcendental equation for the lattice parameter of the equivalent crystal (see 
Ref. 24 for a detailed derivation of this equation) 


NR[e 


-at? i 


+ MR 2 e 


a + ; 


X p -(a + S(r t ))r k 
r k e 
k 


(3) 


where N and M are the ideal number of nearest- neighbors {N) and next-nearest neighbors (M) in 
the equivalent (ideal) crystal, Rj and R 2 are the nearest-neighbor and next-nearest-neighbor dis- 
tances, respectively, in the equivalent crystal, and where p, a and / are ECT parameters that fully 
describe the corresponding atomic species within the context of ECT [24], Eq. 3 is solved for a s , 
the lattice parameter of the equivalent crystal {Rj and R 2 are the nearest neighbor and next nearest 
neighbor distances, respectively, in a crystal with lattice parameter a s ). The right hand side (r.h.s.) 
of this equation, which is computed in terms of the actual distances for the known structure in the 
defect crystal r k , can be interpreted as a measure of the defect (in the sense of how the actual envi- 
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ronment surrounding the reference atom differs from that in its pure crystal) due to the actual 
atomic environment of the reference atom where r k are the distances between these atoms and the 
reference atom. S(r) is a screening function, also to be defined later. The sum runs over all neigh- 
boring atoms (at distances r k ) within a sphere of radius r c , defined as the next-nearest neighbor 
distance in the perfect euqilibrium lattice. 

If this equation admits a solution, that means that the equivalent crystal associated with the 
reference atom is properly defined and its energy can be determined by use of the binding energy 
relationship (demonstrated to accurately represent the UBER) [25,26]: 


= ~EJl + a s )e 


~a s 


(4) 


The scaled lattice parameter a* s is related to the solution of Eq. 3, a s , by means of the fol- 
lowing expression 


* ( a S~ a e) 

a s = c l 1 


(5) 


where a e is the equilibrium lattice parameter, E c is the cohesive energy of the ground state crystal, 
q is the ratio between the equilibrium lattice parameter and the equilibrium Wigner-Seitz radius 
r WSE, an d / is a scaling length to be defined later. Having computed the equivalent crystal, we now 
determine the BFS strain energy contribution of the reference atom as the difference between the 
energy of the atom in the defect crystal and its energy at equilibrium in its ground state crystal: 


e S = E c l-0+rt S )e 


( 6 ) 
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where the subindex S to the scaled lattice parameter denotes the fact that this equivalent crystal is 
related to the strain energy component of this atom in the alloy. The parameters p, a, / and X are 
ECT parameters that fully describe the atomic species in question: p is 2/7-2, where n is the princi- 
pal quantum number, a e is the equilibrium lattice parameter, / is a scaling length and X is a 
screening factor to properly account for screening effects beyond first neighbors (see Ref. 24 for a 
full description of the determination of these parameters). 

The parameter a, which is directly related to the structure of the wave function in the overlap 
region between similar atoms, has been determined in past applications of ECT by requiring that 
the formalism reproduces the experimental value of the single vacancy formation energy [24,25]. 
However, in order to minimize the dependence on experimental or theoretical input for the deter- 
mination of these parameters, we now introduce a novel way of determining the parameter a, 
which provides the additional benefit of determining the range of validity of the approximations 
introduced in writing the perturbation equation for the calculation of the equivalent crystal. As 
mentioned before, the r.h.s. of Eq. 3 can be considered as a ‘measure’ of the defect (in the sense 
that the departure of this quantity from its equilibrium value is representative of the difference 
between the defect crystal and the perfect, equilibrium crystal). For example, the presence of 
vacancies close to the reference atom reduces the value of this quantity, while other defects, like 
interstitials, increase its value. While a solution exists for arbitrarily small values of this quantity, 
this is not so for values much larger than the one corresponding to equilibrium (for which the 
equivalent crystal is the equilibrium crystal itself). There is a maximum value for the r.h.s of Eq. 3 
above which no solution to this equation can be found. We therefore determine a by requiring that 
this maximum defect corresponds to the situation when the strain energy of the reference atom 
equals its cohesive energy. In doing so, a becomes dependent on the other parameters characteriz- 
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ing the atomic species (E c , a e , p, l, X). This redefinition of the ECT parameter a reduces the actual 
necessary input to just the cohesive energy (E c ), the equilibrium lattice parameter (a e ), and the 
bulk modulus (B 0 ) since / is directly related to these three parameters by 


2 E c 

r = 


l2nr W sE E o 


(7) 


where r WSE denotes the equilibrium value of the Wigner-Seitz radius. Furthermore, the screening 
length X [24] is defined as 


TABLE 1 

COMPUTED INPUT PARAMETERS FOR NiAl CALCULATIONS 


LMTO results 

Atom 

Lattice 

Parameter 

a e (A) 

Cohesive 

Energy 

E c (eV) 

Bulk 
Modulus 
B 0 (GPa) 

Vacancy 
Energy 
E v (eV) 

Ni 

2.752 

5.869 

249 

3.0 

Al 

3.190 

3.942 

78 

1.8 

ECT parameters 

Atom 

P 

a 

1 

X 



(A' 1 ) 

(A) 

(A) 

Ni 

6 

3.067 

0.763 

0.2717 

Al 

4 

1.8756 

1.038 

0.3695 

BFS parameters 


^NiAl - -0.0581 A 1 A A1Ni - 

0.0840 A' 1 



LMTO results and ECT parameters for bcc-based Ni and Al, and BFS 
parameters for B2 NiAl. 
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X = 2.81/. 


( 8 ) 


Finally, the screening function S(r) is defined by S(r )= 0 if r<r j, S(r) = (1/2X)(1 - cos \\) if 
r 2 <r<r i (where (3 = n(r-r 1 )l(r 2 -r 1 )) and S(r)= 1/X if r>r 2 , where r j and r 2 are the equilibrium 
nearest- and next-nearest-neighbor distances in the equilibrium crystal [24,25]. 

Chemical energy 

The calculation of the chemical energy involves an environment related to the ground state of 
the reference atom (i.e., the reference atom, as well as its neighbors, are located in equilibrium lat- 
tice sites corresponding to the ground state crystal of the reference atom) although the actual 
chemical distribution of the atoms is conserved. Given the actual chemical environment of a given 
atom (those atoms within next-nearest neighbor distance from it), it is straightforward to deter- 
mine an equivalent chemical environment in a perfect crystal. A set of transformations can be 
defined which assigns, to each site in the perfect crystal, the appropriate chemical species. The 
resulting distribution is then chemically equivalent to that found in the alloy, hi essence, the same 
considerations made for the strain energy apply here, with the only difference being that the 
changes of electron density in the vicinity of the reference atom are now due to changes in the 
atomic identity of the neighboring atoms rather than changes in the atomic locations, hi this sense, 
the same concepts (the existence of an equivalent crystal) apply, only now the reference atom can 
have different energies than that allowed by its own UBER. As mentioned before, a is the only 
parameter within our formalism that carries all the information regarding the electron density in 
the overlap region between a pair of atoms. For dissimilar atoms it is only natural then to ‘perturb’ 
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such an environment by allowing the wave function in that region to be parameterized by a 
slightly different parameter than that used hi monatomic calculations. We therefore define the 
BFS parameter A ab (and A ba ), given by 

a AB = a A + A ba (9) 

when the reference atom is of species A and its neighbor of species B and 

a BA = a B + A ab (10) 

when B is the reference atom and its neighbor is an A atom. The BFS parameter A ab (or A ba ) 
perturbs the pure element parameter a, indicating the mixed nature of the bond. 

Following the guidelines for writing Eq. 3 and taking into account that the neighboring atoms 
are by definition located in equilibrium lattice sites of the reference atom, the trascendental equa- 
tion to be solved for the (chemical) equivalent crystal is simply 

-<«4 + Au hi , 

e + 

k k 

where r j A denotes the nearest neighbor distance in an ideal crystal A and r 2A the corresponding 
next-nearest neighbor distance. All other terms have their previous meanings. The calculation of 
the chemical energy proceeds from solving Eq. 1 1 for a c , the lattice parameter of the equivalent 
crystal for which R 1 and R 2 are the nearest- and next-nearest neighbor distances, respectively,, and 
then determining the scaled lattice parameter a* c given by 


Z 


Pa 
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( 12 ) 


* (a c -a A ) 
a c = ( 1 1— 

1 A 

so that the chemical energy is 

e c = n 1 E c {^ -(! + a c) e *‘1 (13) 

where y= 1 if a* c > 0 and y = -1 if a* c <0 [18], 

Finally, a safeguard should be introduced in the definition of the chemical energy contribution 
to account for those situations where the reference atom does not have full coordination (i.e., the 
number of nearest-neighbors is less than that found in the perfect, equilibrium crystal). Consider 
the reference atom occupying a surface site of an alloy, hi computing the chemical energy, infor- 
mation about the existence of the surface is introduced by the fact that there are not enough atoms 
in the vicinity of the reference atom to correspond to the assumption in Eq. 1 1 . The chemical 
energy obtained in this way would carry information not only on the chemical effect but also on 
the structural effect due to the absence of some neighbors, hi order to completely free the chemi- 
cal energy from this structural information, we reference the previously defined chemical energy 
to a similar structural state but where all the atoms surrounding the reference atom are of the same 
identity as the reference atom. A reference chemical energy e C(S is thus computed in this manner, 
so that the total chemical energy contribution from the reference atom is 

8 c = e c“ e c 0 • (14) 

where e c is computed using Eq. 1 1 but setting the BFS parameters to zero. The energy is then 
computed using Eqs. 12 and 13. 
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Coupling the strain and chemical energy 

As defined, the chemical energy does not depend on the actual geometry of the alloy which is, 
obviously, unrealistic. Therefore the BFS strain and chemical energy terms need to be recoupled 
by properly accounting for the influence of structural effects on the chemical contribution. In 
order to properly describe the asymptotic behavior of such a quantity ( e c should vanish at large 
separations and should drastically increase at small separations), we introduce a coupling function 
g linking both terms. The coupling function introduces this asymptotic dependence by means of 
a* s , the strain scaled lattice parameter (Eq. 5), which can be understood as a measure of the strain 
effect. If a* s = 0 it means that the reference atom finds itself in an environment that resembles 
equilibrium; a positive value of a* s results from average expansions with respect to equilibrium 
and negative values of a* s relate to average compressions. We therefore define the coupling func- 
tion g as: 

g = e~° s (15) 

so that defects that involve expansions reduce the effect of the chemical energy on the total energy 
of formation and viceversa. Fig. 6 displays the effect of the coupling function g on the chemical 
energy, which otherwise would amount to a constant contribution to AH, independent of the size 
of the cell. Summarizing these concepts, the BFS contribution e from an atom i to the energy of 
formation AH is given by 

e = e s + g(e c -e Co ) = e s + ge c (16) 

A graphical representation of the calculation of e is shown in Fig. 7. 
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Fig. 6: Schematic respresenta- 
tion of the strain and chemical 
energy contributions of an atom 
A in an alloy as a function of 
the (strain) equivalent lattice 
parameter a s . The dashed line 
indicates the constant value of 
the chemical energy, which 
amounts to the total contribution 
to the energy of formation only 
when the lattice parameter of the 
alloy coincides with that of the 
equilibrium crystal A. The total 
chemical energy contribution, 
gee, vanishes at large inter- 
atomic separations. For exam- 
ple, if the equivalent (strain) 
lattice parameter is a 1 , we show 
the corresponding strain (e s ) and 
chemical (gee) energy contribu- 
tions. If the equivalent (strain) 
lattice parameter is a e , then e s = 
0 and e c corresponds to its con- 
stant value ( g = 1). 


Applying BFS to the calculation of the energy of formation of a bcc binary alloy 


To demonstrate the ideas presented thus far, we apply the BFS formalism to the calculation of 
the energy of formation of a binary A-B alloy in the B2 structure, as a function of the alloy lattice 
parameter a 0 . Fig. 8 represents the alloy formation process where two pure bcc crystals A and B, 
each with equilibrium lattice parameters a A and a B are reassembled to form a B2 alloy with lattice 
parameter a 0 . hi the abscence of defects, this is a particularly simple example due to its high sym- 
metry and uniformity. All A atoms are equivalent, and so are the B atoms. It is only necessary to 
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compute the energy contribution of one generic A atom and one B atom. In what follows, we will 


compute e A only, since e B would be computed in the same maimer. The different contributions to 
e A are shown in Fig. 9. 

In writing the equation for the calculation of the BFS strain energy (see Eq. 3), we note that in 
the actual alloy an A atom has 8 nearest neighbor B atoms at a distance r j and 6 next-nearest 
neighbor A atoms at a distance a 0 . Following the method for evaluating the strain energy environ- 
ment, we must remember to ‘flip’ the identity of these neighbors to that of the reference atom. The 
BFS strain energy equation for atom A is then 



Fig. 7: Schematic representation of the BFS contributions to the total energy of fonnation. The left hand 
side represents the reference atom (denoted by an arrow) in an alloy. The different terms on the right hand 
side indicate the strain energy (atoms in their actual positions but of the same atomic species as the refer- 
ence atom), the chemical energy temi (atoms in ideal lattice sites) and the reference chemical energy 
(same as before, but with the atoms retaining the original identity of the reference atoms). 
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This equation is trivially solved, with the lattice parameter of the equivalent crystal being just 
a s = a 0 • The corresponding BFS strain energy contribution is 


A j-,A( . , .. , A*. ~a s 

e S = E c 1+ ( 1+ «S ) e 


(18) 


where the scaled lattice parameter of the (strain) equivalent crystal is given by 


A* 

a S 


= Q- 


(«o -a A ) 

T a 


(19) 


A similar calculation is carried out for atom B, replacing p A , l A , a A and X A with p B , l B . a B and 
k B in Eqs. 17-19. From the calculation of the strain energy contributions for atoms A and B one 
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Fig. 8: Crystals A and B (with lattice parameter a A and a B ) are 
reassembled to form a B2 A-B alloy with lattice parameter a 0 . 
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Fig.9: Schematic representation of the contribution of an A atom (center of the cube) to the energy of 
formation of B2 A-B compound with lattice parameter a 0 . The left hand side represents the actual B2 
structure. The first term in the left hand side indicates the strain energy environment (all atoms are of 
type A, in the lattice of the alloy and the second and third term (between brackets) indicate the chemical 
energy environment and the reference chemical energy environment, respectively. 


would also derive the corresponding coupling functions g A and g B , according to Eq. 15. 

Fig. 10 shows the environment assumed in the calculation of the chemical energy. Since in the 
example given, both pure crystals are of the same crystallographic structure as the alloy, the cal- 
culation is particularly simple, and it only requires noticing that the lattice parameter for each case 
is the equilibrium lattice parameter of the pure crystals (a A and a B ). 

The BFS equation for the calculation of the chemical energy contribution for atom A is 



0 
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5 ? 


o 

■o 
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A {a A ) 


B (a B ) 


Fig. 10: Environment seen by atoms A and B, respectively, for the calculation of the BFS chemical energy 
contribution. Note that the actual atomic composition of the alloy is preserved, but the lattice parameters 
used for each case corresponds to that of the equilibrium pure crystals A and B, respectively. 
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where r A is the equilibrium nearest neighbor distance in a crystal A. Note that atom A interacts 
with 8 nearest neighbors of species B (thus the factor a A +A BA in the first term of the r.h.s. of Eq. 
20) and that all 6 of its next-nearest neighbors are of the same atomic species A. This equation is 
solved for the equivalent lattice parameter a c and the chemical energy contribution is 


e c = 


y a4 


,, a* ■ 

(1 + a c )e 


( 21 ) 


where the scaled lattice parameter of the (chemical) equivalent crystal is given by 


A* 

a C 


, A 

{a c 


«a) 


( 22 ) 


and where y A = 1 or -1, depending on whether the scaled lattice parameter is positive or negative. 
A similar calculation is performed for atom B, replacing the appropriate parameters. In particular, 
note that the interaction parameter in the exponential in the first term of the r.h.s. of Eq. 20 is now 
0 Cg+A AB , as all the nearest neighbors of the reference B atom are of species A (Fig. 10). 

In this particularly simple example, it is not necessary to compute the reference chemical 
energy (<? c ) for either type of atom, as the corresponding environments correspond to perfect 
equilibrium A and B crystals already. Therefore, the total chemical energy contributions are 

A A j B B 

£ c = e c ailQ £ c = e c . 
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We finally write the following expression for the energy of formation AH of the B2 structure 

A-B: 

\ jt 1 / -A -A B B+ ,-av 

AH = ~(e s + g A e c + e s + g B e c ) (23) 

While this expression applies to this stoichiometric binary system characterized by one atom 
of each atomic species, in more complex situations the calculation might involve considering 
more non-equivalent atoms. In this example the Ni atoms locate themselves in one simple cubic 
sublattice and the A1 atoms in the other. As such, all Ni or A1 atoms are respectively identical and 
contribute the same amount to the energy of formation. Therefore, the BFS calculations have to 
be performed for just one atom of each species. The simplicity of this calculation allows for a 
straighforward determination of the BFS parameters hi terms of the energy of formation and equi- 
librium lattice parameter of the binary base alloy (NiAl). In the following subsection, we will 
refer to the scheme used in this work for parameterizing the binary system at hand. 

Parameterization of the BFS Method 

The previous example indicates that beyond the approximate scheme developed to simulate 
the process of alloy formation (i.e. separating strain and chemical contributions), the power and 
accuracy of the method are heavily dependent on the parameters used both for the individual ele- 
ments (p, I, < 7 e , E c , X) and those used for each type of binary combination that might appeal* in the 
system (A ab and A ba ). Also, this example highlights the extreme simplicity of the method when 
the crystallographic structure of the alloy is the same as that of the constituents (in the example, 
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both A and B were bcc elements and the alloy formed is also bcc-based). While in general it is not 
known what is the phase structure of an alloy for a certain composition, some basic knowledge 
can be used a priori to simplify the calculations. 

Take, for example, the system to be studied hi this work: (3-NiAl. This binary alloy is known 
to have the B2 structure, in spite of the fact that each of its constituents are fee elements. Assum- 
ing the existence of a bcc phase of Ni as well as for A1 would reduce the complexity of the calcu- 
lation to a level similar to the example just discussed. If the BFS method relied only on 
experimental input for the determination of the necessary parameters listed above, it would make 
this approach impossible to implement under most circumstances. However, first-principles 
methods can be used in a straightforward maimer to determine the physical properties of these 
elements ( bcc-Ni and bcc-Al) even if they do not exist naturally. In doing so, not only do we 
eliminate the dependence of the BFS method on sometimes uncertain or non-existent experimen- 
tal input, but we also generate a consistent approach. This approach simplifies the numerical 
complexity of the problem and also introduces a systematic procedure for the generation of the 
necessary parameterization, since they would all be computed by the same method, with the same 
level of approximation. 

To determine these parameters, we used the Lmear-Muffin-Tin-Orbital (LMTO) [27] method 
in the Atomic Sphere Approximation (ASA). This scheme was used to calculate the equilibrium 
properties ( cohesive energy E c , equilibrium lattice parameter a e , and bulk modulus B 0 ) of the ele- 
mental solids in the same crystal symmetry as that of the compound to be studied (bcc). This set 
of parameters is accurately described by the Local Density Approximation [28]. In this work, we 
generated parameters for Ni, Al, Ti, Cr and Cu in the bcc symmetry, as well as the B2 structures 
of each of the pairs formed (only B2 NiAl and NiTi exist in nature) for the determination of the 
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BFS parameters A ab and A ba . The parameters for Ti, Cr and Cu will be used and validated in 
subsequent papers in this series. 

The LMTO method uses a minimal basis set. For this work, we only used 5, p and d orbitals. 
All calculations were done with equivalent sampling of the Brillouin zone using, for the bcc lat- 
tice, 120 ^-points in the irreducible wedge. The basis set and sampling of the Brillouin zone used 
for the calculation of the B2 compounds was equivalent to that used for the pure elements. We 
refer the reader to Refs. 27-28 for a detailed discussion of the LMTO method. The pure element 
LMTO-generated parameters p, l, E c , a e and the resulting value of a are displayed in Table 1 for 
Ni and Al. 

The simplicity of the BFS formalism allows for a straighforward method for determining the 
BFS parameters A ab and A ba by an analytical procedure that not only eliminates numerical 
errors but also sheds light on the range of validity of the resulting parameters. Leaving the details 
of the determination of the BFS parameters to Appendix 1, we list in Table 2 the resulting values 
for A NiA1 and A A1Ni . 

Other theoretical methods applied to intermetallics 

As mentioned in the Introduction, there has been substantial progress in the last two decades 
in the area of theoretical techniques and their relevance in dealing with systems of practical 
importance. Current techniques range from those that impose an empirical parameterization of the 
constants in Ising Hamiltonians [29-33] or those that parameterize experimental databases [34], 
families of semiempirical methods based on the determination of interatomic potentials [11-15], 
to all-electron first principles calculations [35-44], Recent first-principles calculations include the 
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full-potential linearized augmented Slater-type orbital method (LASTO) used by Alatalo et al., in 
their analysis of the complex Zr-Al phase diagram [40], the first-principles discrete variational 
method study of disordered Fe-V alloys by Krause et al. [41], a detailed study of composition and 
phase diagrams of Cu-Au, Ag-Au, Cu-Ag and Ni-Au by Ozolins et al. [42], and the determination 
of formation energies for Fe 3 Al by Mayer et al. [43]. These, and other similar studies, concen- 
trate mostly on monatomic [12,44] or binary systems. 

In terms of semiempirical techniques, current work on NiAl alloys include the extensive anal- 
ysis of this system by Farkas et al., in a series of recent articles where the embedded atom method 
is modified for its use in this system [45-46]. Other techniques currently used include Finnis-Sin- 
clair potentials specifically determined for NiAl [47], tight-binding methods [12,44] (still mostly 
restricted to single element systems), and potentials determined from first-principles calculations. 
This last group includes the work of Zhang et al. in Fe-Al and Ni-Al alloys [48], and the study of 
transition metal aluminides by Moriarty and Widom [49]. 

Other relevant work in this area includes the extension of Miedema’s semiempirical model for 
point defects in B2 compounds by Bakker et al. [34] and the development of semiempirical N- 
body non-central potentials for FeAl alloys by Besson and Morillo [50]. 

In terms of ternary NiAl+X systems, very little progress has been made relative to that made 
for binary alloys, although it is expected that the development of accurate potentials and efficient 
computational methods, will generate growth in this area in the near future. Most notably, recent 
theoretical work of the influence of ternary additions to NiAl includes the Monte Carlo analysis of 
NiAl+Ti alloys by Sumin et al. [51], and the statistical theory of ordering combined with the elec- 
tronic theory of alloys in the pseudopotential approximation by Mekhrasov et al. [52], 

As we increase the complexity of the system at hand, as will clearly be the case in the subse- 
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quent papers in this series, we will have little if any previous theoretical work for comparison. 
Once we move beyond ternary systems, all our results will become purely of a predictive nature, 
relying on experimental verification. 

RESULTS AND DISCUSSION 

The BFS calculations for B2 NiAl were performed on a 72-atom cell, allowing for atomic 
relaxation induced by the presence of vacancies and substitutional atoms. However, relaxations 
were allowed only at a global level in the sense that the lattice parameter of the cell as a whole 
was varied until the energy was minimized. A more detailed study, which we consider unneces- 
sary for the goals of this work, would include the possibility of local atomic relaxations in the 
vicinity of defects. This issue, and its consequences on the phase structure, will be dealt with in 
forthcoming papers. 

In a cell of this size there is a large number of possible distributions of Ni and A1 atoms for 
any given alloy composition. Knowing the behavior of each and every one of them would be 
unnecessary (in this case, considering all the possible distributions of Ni and A1 atoms as well as 
vacancies, in a 72-site cell), as most of the information needed for determining the defect struc- 
ture can be obtained from a few selected situations where the essential defect structures are con- 
sidered. Therefore, instead of searching for an absolute energy minimum for a given number of Ni 
and A1 atoms located in every available lattice site with any specific short and long order scheme, 
we chose to construct a selected, but sufficiently large number of 'candidate’ distributions to 
obtain information on the energetics of the system close to the ground state. These configurations 
are still, to a certain extent, arbitrary, and may not necessarily include the ground state for a given 
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composition. However, it is clear that when creating a catalogue of selected configurations which 
consider the most likely atomic distributions, one would obtain not only the necessary informa- 
tion (i.e., the lowest energy state for each concentration) but also additional information on gen- 
eral trends in defect structures. While most of the configurations included in the catalogue are 
experimentally unaccessible, they serve the purpose of determining trends in the bulk properties 
of the actual alloys. Moreover, if the set of configurations sampled is sufficiently large and the 
structures are properly chosen, respecting the symmetries that characterize the system, one would 
expect to find the ground state or states sufficiently close to it at each composition. 

As a consequence, an additional advantage of this approach is that it enables us to determine 
an ‘energy spectrum’, showing how alternative ordering schemes compare energetically with the 
lowest energy state for each composition. The magnitude of the energy gaps between alternative 
distributions and the ground state should be proportional to the probability of finding these config- 
urations in the actual alloy. Figure 1 1 includes representative samples of the catalogue of configu- 
rations used in this work, first showing the basic 72-atom cell (Fig. 1 1 .a) followed by a few cases 
that include substitutional defects, antistructure atoms and vacancies (Fig. ll.b-d). 

Based on stoichiometric NiAl, Appendix 2 introduces a simple way of characterizing each one 
of the configurations included in the catalogue of structures we examined. Due to the large num- 
ber of configurations included, only a small but representative number of them are graphically 
shown in Fig. 11. The rest are described in a notational shorthand discussed in Appendix 2. 

The analytical procedure is then straightforward. The BFS method is applied to each one of 
the configurations included in the catalogue in order to determine its energy of formation as well 
as the equilibrium lattice parameter (for which the total energy of formation is minimized). 
Appendix 2 also includes these BFS results for each configuration. 
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Fig. 11: (a) The computational cell used to generate the catalogue shown in Appendix 2. (b) Sam- 
ple of a NiAl cell showing two site exchanges (see [4] in the catalogue). 
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Fig. 11 (continued): (c) Sample of a Ni 52 . 7 sAl 47.12 alloy where A1 atoms at lattice sites 38, 41 and 
46 are replaced by Ni atoms (see [25] in the catalogue in Appendix 2). (d) Sample of a 
Ni 47 . 83 Al 52.17 alloy, where the A1 atom at lattice site 42 is replaced by a vacancy, the Ni atom at 
lattice site 30 is replaced by an A1 atom, and the Ni atom at lattice site 31 is replaced by a vacancy 
(see [38] in the catalogue). 
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Fig. 12 displays some of the results for representative configurations in the catalogue. The 
results are shown in the form of an energy spectrum, which shows a clear comparison of the ener- 
getics of the various atomic arrangements. Each column hi this spectrum corresponds to specific 
compositions of Ni 5 o +x Al 5 o_ x alloys, denoted by the Ni atomic concentration. Each energy level is 
identified by the catalogue number used in Appendix 2. The lowest energy configuration (lowest 
level in the first column in Fig. 12) corresponds to the perfect B2 structure, where all Ni and all A1 
atoms populate their own interpenetrating sublattices. Higher energy configurations represent a 
wide variety of defects. Moving up hi the spectrum, we first find a group of levels very close in 
energy which correspond to the presence of antistructure atoms, where a pair of Ni and A1 atoms 
switch places. The small differences in energy between these levels relates to the relative position 
of the two antistructure atoms (being smallest when they are nearest-neighbors). After another 
comparable gap in energy, we found states with two antistructure pans and so on. Similar situa- 
tions are found at higher Ni concentrations. The lowest energy configuration for a given concen- 
tration in Ni-rich alloys always corresponds to a pure substitutional state and higher energy states 
include additional swaps of atoms between Ni and A1 sites. For the highest Ni concentration 
shown (55.55 at. % Ni), the purely substitutional case distinguishes between the configuration 
where the substitutional atoms are as far away from each other as possible, to configurations with 
slightly higher energy where the substitutional atoms display a clustering tendency. 

This spectrum emphasizes both the completeness of the catalogue selected (i.e. the lowest 
energy state is properly identified for every concentration) as well as its shortcomings: some high 
energy states (2 for 50 at. % Ni; 19 for 51.39 at. % Ni), characterized by clustering of antistructure 
atoms) are missing at higher concentrations. On the other hand, the high energy state 27 for 55.55 
at. % Ni is comparable to 8 for 50 at. % Ni, both including specific alignment of the antistructure 


NAS A/TM— 1998-208820 


38 




50.00 51.39 52.78 54.17 55.55 

at. % Ni 


Fig. 12: ‘Energy spectrum’ for substitutional defect structures in Ni-rich NiAl alloys. Each column 
represents a specific binary alloy concentration. The levels represent the energy of formation of con- 
figurations specified in the catalogue (the numbers by each level identify the configuration in Appen- 
dix 2). All the configurations included in this Figure correspond to substitutional alloys exclusively 
(no vacant sites). In all cases, the lowest energy configuration corresponds to Ni atoms occupying 
available A1 sites. Higher energy states include various antistructure defects. 
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atoms in specific planes and directions (see Appendix 2). The smaller energy gap between this 
type of configuration and the lowest energy state for higher Ni concentrations hints at the possibil- 
ity that this type of defect might occur in the actual alloy. 

Fig. 13 shows similar results for Al-rich alloys, indicating not only that the lowest energy 
states arise from the presence of Ni vacancies but also that selected antistructure-vacancy combi- 
nations could also be present. Experimental evidence confirms these results [4-9]. 

The determination of the lattice parameter for each configuration is done by minimizing the 
energy of the cell with respect to the lattice spacing of the cell as a whole but not allowing for 
individual anisotropic relaxations. While performing detailed calculations including individual 
atomic relaxations is well within the computational efficiency of the method, such lengthy calcu- 
lations would become impractical for higher order alloys. We therefore compromise between 
accuracy and simplicity assuming that whole-cell isotropic relaxations will be sufficient for 
describing the energetics of these alloys. 

The wealth of information embedded in Fig. 12 can be used to gain further insight in the 
dependence of the lattice parameter a(x) of Ni 5 Q +x Al 50 _ x alloys. While it is natural to determine 
a(x) from the values of a corresponding to the lowest energy states for each concentration, as 
depicted in Fig. 12, it is interesting to first visualize the lattice parameter of alternative, higher 
energy, structures by choosing characteristic values for each ‘band’ (group of energy levels close 
in energy) in Fig. 12. 

Fig. 14 displays the lattice parameter as a function of Ni concentration for a portion of the set 
of 'candidate’ configurations, regardless of the energy of the configuration. Each point corre- 
sponds to a distinct arrangement of Ni and A1 atoms as well as vacant Ni or A1 sites in the calcula- 
tional cell. We focus our attention on three particular regimes that can be seen in Fig. 14: an 
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inverted V-shaped set of filled symbols with its apex at the stoichiometric composition, a line of 
circles indicating structures with increasing lattice parameter with increasing A1 content, and a 
large number of configurations enclosed by the boundary line consisting of filled symbols. The 
configurations denoted by filled symbols represent the lowest energy configurations at each com- 
position. The ascending line of circles corresponds to configurations where no vacancies are 
present and increasing numbers of A1 atoms occupy Ni sites, a typical substitutional defect struc- 
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49.29 48.57 47.83 47.06 46.27 

at. % Ni 

Fig. 13 : Energy spectrum for Al-rich NiAl alloys containing vacancies. For x Ni = 49.29, the low- 
est energy state corresponds to an A1 atom in a Ni site and a nearest-neighbor vacant A1 site (34), 
a vacant Ni site (33) and an A1 atom in a Ni site and a distant vacant A1 site (35), in order of 
increasing energy. For x Ni = 48.57, the lowest (two-vacancy) state corresponds to an A1 anti- 
structure atom and a nearest-neighbor pah of vacancies in Ni and A1 sites (38) with almost the 
same energy of formation as a configuration with two vacant Ni sites (37). For higher A1 concen- 
trations, the lowest energy state is always one with vacant Ni sites (40, 43, 46). 
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Types of Defect Structures 
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(3-NiAl 
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□ 

Ni -antisite 

X 

Ni and A1 
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A 

Ni vacancies 

<> 

A1 vacancies 
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Ni and A1 
vacancies 

o 

Ni -antisite 
and Ni vacancies 

Al-antisite and 
A1 vacancies 


o 

Fig. 14: Lattice parameter (in A) as a function of Ni concentration for various candidate defect 
structures. Different symbols indicate alloys with different types of general defect structures. Sim- 
ilar symbols for the same composition correspond to different distributions of the same defect 
within the computational cell. The shaded symbols correspond to those configurations with the 
lowest energy for that composition. 
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ture. While substitutional states have the lowest energy for each concentration in the Ni-rich side, 
the opposite is true for the 45-50 at. % Ni concentration range. The cluster of points enclosed by 
the inverted V-shaped boundary correspond to configurations characterized by random distribu- 
tions of vacancies in A1 or Ni sites as well as random exchanges of Ni and A1 atoms from their 
equilibrium stoichiometric locations. Stoichiometric [1- NiAl is noted by an asterisk. 

A comparison of the modeled results with experimental data for the lattice parameter of vari- 
ous NiAl alloys is displayed in Fig. 15. The theoretical results, shown as a continuous line, consist 
of the calculated lattice parameter and density for the lowest energy configuration observed at 
each composition. The lattice parameter results from the lowest energy configurations on the Ni- 
rich side of stoichiometry, show excellent agreement with experimental data (Fig. 14). The linear 
regime in the density (Fig. 15. a) and lattice parameter (Fig. 15.b) vs. Ni concentration for the 
range 50-60 at.% Ni is almost exactly reproduced by our calculations. A recent review paper by 
Noebe et al. [53] provides a convenient linear description of the available experimental results, 
both for the lattice parameter a and density p vs. Ni concentration. To highlight the agreement 
between theory and experiment, we therefore compare those expressions with ones corresponding 
to the BFS results. For comparison purposes, both quantities were normalized to the stoichiomet- 
ric NiAl values {a 0 and p 0 , for the lattice parameter and density, respectively), as there is a slight 

O O 

difference of 1.2 % between the experimental {a 0 = 2.887 A) and LMTO values {cij heor = 2.85 A) 
for the lattice parameter of the stoichiometric B2 NiAl structure. The expressions for experimen- 
tal and theoretical results on the Ni-rich side of stoichiometry are 

f-2-1 = 0.5339 + 0.009322.% (24) 

^Po 'Exp 

and 
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(a/a NiA1 ) (p/pNiAl) 




Fig. 15 : (a) Density and (b) 
lattice parameter of binary 
alloys as a function of Ni con- 
centration, normalized to 
their stoichiometric values. 
The solid symbols denote 
experimental results from dif- 
ferent investigators (see Ref. 
53). The lines denote the BFS 
predicted properties (see 
text), corresponding to the 
structure with the lowest 
energy state for each compo- 
sition. 


at. % Ni 
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£ ) 

PoJ 


= 0.5324 + 0.009352%. 


Theor 


for the density (relative to the stoichiometric NiAl density) and 


(25) 


<7 A 

°0 'Exp 


1.0384- 0.00079 14%,. 


(26) 


and 

f— ") = 1.0396-0.0007689%, (27) 

^ a 0 'Theor 

for the lattice parameter (relative to the stoichiometric NiAl value) for Ni concentrations, x Ni , 
between 50 and 60 at. % . 

On the Al-rich side of stoichiometry a more complicated defect structure is observed. The 
lowest energy states on the Al-rich side correspond to arrangements of vacancies with no change 
in site occupancy for the Ni and A1 atoms with respect to their original sublattice in NiAl. In the 
lowest energy configurations, the Ni vacant sites are at next-nearest neighbor distance from each 
other. In other words, clustering of the defects results in the lowest energies. Once again, the ensu- 
ing values for the lattice parameter and the density as a function of Ni concentration (closed trian- 
gles in Fig. 14) can be described by linear relationships and compared to the experimental values 
[19]: 

(£) = 0.02 Sfli (28) 

y V0'Exp 

and 
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( 29 ) 


(30) 


— = 0.840261 + 0.003 184jc w (31) 

^ Cl 0 'Theor 

for the lattice parameter (relative to the stoichiometric NiAl equilibrium value) for x Nj between 45 
and 50 at. %. 

Wliile not directly apparent from Fig. 14, this work also showed that locally ordered arrange- 
ments of defects were energetically preferred to random distributions of point defects. For exam- 
ple, in Fig. 14 at 47.82 at. % Ni there are an open and filled triangle of very similar lattice 
parameter. The difference in these points is that the input structure of the lowest energy configura- 
tion (solid triangle) was composed of a greater degree of ordering of Ni-vacancies. While overall 
structures composed of Ni-vacancies had the lowest energy at a given composition, a similar trend 
of lower energy with increasing ordering of point defects was observed within the other types of 
defect arrangements as well. For example, within the family of configurations composed of Al- 
vacancies plus Ni-vacancies, the configurations (hexagons) with lattice parameters closest to 
those of the ground state consisted of locally ordered arrangements of vacancies (i.e. an A1 
vacancy shared by two Ni-vacancies at the nearest neighbor distance) while those configurations 
with higher energy consisted of more random distributions of point defects. 

This also hints at the possibility that the actual defect arrangements on the Al-rich side of sto- 
ichiometry may be more complicated than normally assumed. A complete treatment of tempera- 
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ture and local relaxation effects would be necessary to absolutely rale out the existence of these 
more complicated structures consisting of ordered arrangements of Ni and A1 vacancies in Al-rich 
alloys. But even that would not rale out the possibility that these structures may actually be 
observed as metastable states since their energy is not that different from defect structures com- 
posed solely of Ni-vacancies. 

CONCLUSIONS 

We introduced a formalism that combines the type of results achieved by first-principles 
methods with the convenience of setniempirical methods for a simple and straightforward analy- 
sis of situations that cannot be treated in a similar fashion by either individual technique. An intu- 
itive as well as a formal description of the BFS method was presented, with the goal of 
familiarizing the reader with the concepts to be used in this and subsequent papers in this series on 
NiAl-based alloy design. Beyond the operational equations, we introduced other elements that 
will be heavily used when dealing with higher order alloys. These include the definition of a cata- 
logue of predetermined configurations (as opposed to lengthy computational algorithms for the 
search of minimum energy structures), the use of energy spectra to identify trends and to gain 
insight on metastable structures (as opposed to just dealing with lowest-energy states), and the 
description of the observed behavior in terms of BFS strain and chemical energies. 

As our first test of the BFS method we analyzed the defect structure of binary NiAl alloys. 
Our theoretical modeling of the defect structure of [3-NiAl alloys is fully consistent with the con- 
clusions drawn from experimental evidence indicating that substitutional defects dominate in Mi- 
nch alloys and that the presence of vacancies is responsible for the sharp decrease in lattice 
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parameter with increasing A1 content in Al-rich alloys. Moreover, the results hint to the possibil- 
ity that a more complex defect structure exists: the presence of vacancies both in the Ni and A1 
sublattice in Al-rich alloys and a preference for clustering or local ordering instead of a random 
distribution of point defects with any deviation in stoichiometry. In future papers we will build on 
these results by modeling ternary and higher order alloys. 
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Appendix A - BFS Parameterization 


An advantage of the BFS method is that, due to its simple formulation, it allows for a straight- 
forward (analytical) determination of the BFS A ab and A ba interaction parameters, therefore 
avoiding numerical uncertainties inherent to any numerical fitting procedure. Moreover, the input 
data used, whether it is obtained experimentally or from other theoretical calculations, ‘localizes’ 
the accuracy of the ensuing BFS predictions for those alloys in the vicinity of the phase diagram 
of the ordered structure used as input. 

For simplicity, we reduce the following derivation to the case where the ordered structure cor- 
responds to a cubic lattice characterized by a single lattice parameter a 0 (simple bcc or fee alloys 
with no tetragonal distortion). 

Consider an alloy A-B, where due to the symmetry of the structure, there are N x non-equiva- 
lent atoms of species X (X = A, B) and n x denotes the multiplicity of the /th non-equivalent atom 

of species X, so that 

N x 

I X "x, = N c (Al) 

x I. 1 

where N c is the total number of atoms in the cell, hi this context, two atoms are deemed to be 
equivalent if they have the same environment within a sphere of radius R, the nearest-neighbor 
distance in the ground state crystal. In this case, the two conditions used to determine the BFS 
parameters consist of exactly reproducing the heat of formation A H 0 of the ordered structure and 
the corresponding lattice parameter, which are determined via LMTO calculations for a given 
structure [27,28]. In this work, we used the B2 NiAl base alloy as the basis for the LMTO calcula- 
tion. 
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The two conditions are 


N, 


1 I j^«*, = Atf 0 

X t ,,"c 

and 


(A2) 


N 


X n X^ e X: 


'ST' 'ST' _A ' 

V N r da 

x i = 1 


= 0 


(A3) 


where represents the BFS contribution to the energy of formation, given by 


e Xi ~ e Xi + 8x( e x~ e x) 

In Eq. A4, e s x is the BFS strain energy and z c x = e c x -e x is the BFS chemical energy, where e x is 
the reference energy. 

For bulk ordered alloys, such as B2 NiAl, the BFS strain energy is the same for all atoms of 
the same species and it is uniquely determined by the input value of the lattice parameter a 0 , 


e x, ~ Eel 1 _ 1 + a Xi )£x ; j 
where 


(A5) 


a S x. = P x( a O~ a e) ( A6 ) 

where E* and a x e are the cohesive energy and equilibrium lattice parameter, respectively, for 
atoms of species X and \i x j 1 , where q is a structure constant ( q ^ for bcc elements) and i x 
is a scaling length for species X [24,25], The coupling function included in Eq. A4 is given by 

8x, = ( A7 ) 


The BFS energy depends on the lattice parameter of the alloy structure only via the BFS 
strain energy and glue, therefore, Eqs. A2-A3 can be written as 
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X X aj I e x l ( a o) + <?x/ a o) e x, 1 _ ^// 0 

X i= 1*1 J 


(A8) 


and 


/v 


v V 


+ e 






'X, 


<i<7 


= 0. 


(A9) 


Eq. A8 and A9 can then be written as 


n y 


n 


Wv 


XX w £x,( a o) e x, _ X /v e x( a o) 

x i= 1 c 

and 


c 


/7 


X: S 


X /=1^ C 


(A10) 


„ At x 

A, X ; ~ , . c 


E Z W Px£x,( a o) 8 X; _ Z Z W E C$X a Xi( a o)gx( a o) 

Xi = 1 c */=i c 


(All) 

If we concentrate now only on binary alloys that form fee or bcc ordered structures character- 


ized by a single lattice parameter (z.t 2 , l i 0 , B2, B32, etc.), then n A + n B = n c . Simple expressions 
can then be obtained for the BFS chemical energies d, and 4 : 


C ~ ^2 

C — 

«a.?a 0) (Pb-Pa) 

and 

c ^2 ~Pa^i 

g — 

u b 8 b } (P B _ Pa) 

where 


(A12) 


(A13) 


§i = (n A + n B )AH 0 - (n A e S A (a 0 ) + n B e S B (a 0 )) (A14) 

and 

8 2 = Zj^xP x^cSx ( a o) (A15) 

x 

with g x ] = g x (a 0 ) . The BFS chemical energies can then be determined with Eqs. A12 and A13, so 
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that we can then search for the set of parameters ( a ab , a ba ) that simultaneously satisfy these con- 


ditions. This is done by starting with the use of the following expression for the BFS chemical 
energy in terms of the equivalent (chemical) lattice parameter 

4 = Yx4jl-(1 + 4V 4 j (A16) 

where y = l if a x > o and y = -t otherwise. The scaled lattice parameter for the chemical energy 
a c x * , given by 

4* = m%- 4) ( Al7 ) 

is related to the BFS parameters {a ab , a ba ) by means of the BFS equation for the chemical energy 

~( a x + Ax) r x, 

e + 

k k 


r\PX — CCv/?i 

NR\ e 


+ MR P 2 x e 


- a y + r— \R 


= i n 


Px 
Xk r X t 


'Mx/x X 2 e 


*X + Ax 


X x 


(A 18) 


where N(M) is the number of nearest-neighbors (next-nearest-neighbors) in the equivalent crystal 
of species X, /q = ca c x , r 2 = a x , r l = ca x e and r x> = a x e (c = 7 for bcc). If we define 


C»1 f c> 

ca y ( 1 \ ca y 

n l ft '“T* 1 * M ( x 4T 

Qx= Y x ' + ~w\ e \°' + T x e 


X - AT Px 

c h ~ N xx r x l e 


(A19) 

(A20) 


and 


x 

q 2 


~ M xx r x/ 


ax+ x x r*2 


(A21) 


and considering the typical magnitude of the exponent in Eq.A21, it is reasonable to make the 
approximation r x = only in that term, then the l.h.s. of the BFS equation (Eq. A 18) reads 
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X X 


Qx~ <h - <h ~ N X Y r x l e 


Px ~( ax + r X \ 


+ M xy r x e 


Px V k x J 2 


(A22) 


Let 


x 

< Id 


p Y a , p Y 

= N XY r x e + M vv r v e 


(°* + 4> 


XY X, 


(A23) 


so that 


XXX - A rx r x 1 


Qx~<h ~<h = ( id e 


(A24) 


then the BFS parameter a yx is given by 

a _ 1 i \Qx^ c h ~ c l2[ ,,,,, 

ra = 'T | \ (A25) 

This result is exact if second-neighbor contributions are not taken into account, or if for any 
particular structure m xy = o for all X and F, as is the case in this work. The B2 structure is such 
that any given atom has an atom of its own species as a next-nearest-neighbor, i.e. M xx = 6 and 
m yx = o . With the exception of the numerical solution of Eq. (A 1 6), the procedure for the deter- 
mination of a yx using Eq. (A25) is straightforward and simple. Moreover, it can be easily shown 
that in most cases a quadratic approximation to the Rydberg function (l +z)e~ 2 suffices to guaran- 
tee accuracy up to 10 % of the exact results, with the added advantage of a completely analytical 
determination of the BFS parameters a ab and a ba . The parameters used in this work were 
obtained by following the procedure described in this Appendix, including second neighbor inter- 
actions and a numerical solution of the the transcendental equations involved. 
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APPENDIX 2 


In this appendix, we list the catalogue of configurations used in this work for different 
atomic distributions in a 72-atom cell, shown in Fig. 11. The B2 NiAl alloy corresponds to the 
atomic distribution shown in Fig. 1 l.a, where Ni atoms are denoted by black disks (labeled 1, 2, 
3, ...) and the A1 atoms are denoted with open circles (labeled 13, 14, 15, ...). A configuration is 
defined by changing the occupancy of these sites by exchanging atomic positions, substituting 
atoms of one species by atoms of the other species, or by introducing vacancies. The following 
table lists most of the configurations used in this work (only those used to construct the energy 
spectrum shown in Figs. 12 and 13). For example, a Ni atom in site n substituting for an A1 atom 
in site m is denoted with Ni n -> Al m ; V n denotes a vacancy in site n and Al n <-> Ni m denotes an 
exchange of A1 and Ni atoms between sites n and m. The configurations are listed according to 
their Ni concentration. The table also includes the values of the energy of formation of the con- 
figuration as well as its equilibrium lattice parameter (obtained by minimizing the energy of the 
cell with respect to a). 
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x Ni 

Configuration 

AH 

(eV/atom) 

a 

(A) 

50.00 

1. B2 NiAl 

-0.60310 

2.8500 


2 . Ni 3 i <c-> A. 1 ^ 2 * 

-0.50816 

2.8588 


3- Ni 30; 32 <->Al 41)43 

-0.41377 

2.8678 


4- Ni 30; 32 <-> Al 42; 43 

-0.41261 

2.8679 


3- Ni 30 ; 3i <-> Al 42;43 

-0.43544 

2.8656 


6 - Nij _ 9 <-> Al 13 _ 24 

-0.07384 

2.9006 


7- Ni 13 5 7 911 <-> A1i 3) 15 )17)19;21;2 3 

-0.16910 

2.8914 


8 - Alg5_68 <-> Ni 53 _ 5 6 

-0.31271 

2.8773 


9- Ni 30,31 <->Al 18,42 

-0.44184 

2.8650 


10. Ni 3 i <-> Al 42 

-0.50814 

2.8588 


11. Ni 30 <-> Al 46 

-0.49795 

2.8597 


12. Ni 3 o, 31 , 3 ^ 35 , 5 ^ 55 , 53,59 <-> A 1 14 > i 7 j 19 j 22 > 62 > 6 5 56 7 >70 

-0.11404 

2.8964 


13. Ni3 0 ! 3i j 34 j 35 j 54 j 5 5j58j4 9 <->A1 i 3 j i 5 j 21 > 2 3 ; 6 i ; 6 3 ; 6 9 j71 

-0.06528 

2.9011 


14- Ni 50 54>58 -> A1 + Al 62;66>70 ->Ni 

-0.09057 

2.9050 


15. Ni 30 ->A 1 + Ni 2 , 5 , 7 , 10,50,53,55,58 <:> All3,l4,17,18,37,38,41,42 

-0.06414 

2.9020 

51.39 

16. Al 42 -> Ni 

-0.59223 

2.8457 


17. Al 42 -> Ni + [12] 

-0.12280 

2.8902 


18. Al 42 -> Ni + [13] 

-0.06828 

28955 


19. A1i 8 ; 3 o , 42 -> Ni 

-0.50607 

2.8537 

52.78 

20. Al 42 ,43 -> Ni 

-0.57990 

2.8418 


21. A1 42j44 -> Ni 

-0.58163 

2.8416 

54.17 

22. Al 38;42;46 -> Ni 

-0.56609 

2.8380 


23. Al 42 _ 44 -> Ni 

-0.56781 

2.8378 


24. Al 42;46;47 -> Ni 

-0.56759 

2.8379 


25. A1 38 ; 4 i ;46 -> Ni 

-0.53671 

2.8259 
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x Ni 

Configuration 

AH 

(eV/atom) 

a 

(A) 

55.55 

26. Al 38) 4 0;61;63 -> Ni 

-0.56044 

2.8335 


27. Al 65 _ 68 -> Ni + Al 61 _ 64 <-> Ni 53 _ 56 

-0.33123 

2.8541 


28. A1 181 9 ; 42 j43 -> Ni 

-0.55365 

2.8342 


29. Al 38 39 4 1)4 2 -> Ni 

-0.55547 

2.8340 


30. Al 19;2 3 ; 42,46 _> 

-0.55702 

2.8338 


31- Al 18;2 2,23,46 _> Ni 

-0.55523 

2.8340 


32. Al 65 _ 68 -> Ni 

-0.55418 

2.8341 

Al-rich Ni (50 _ x) A 1 (50+x) 

x Ni 

Configuration 

AH 

(eV/atom) 

a 

(A) 

49.29 

33. Ni 31 -> V 

-0.43876 

2.8426 


34. Ni 31 -> A1 + Al 42 ->V 

-0.45149 

2.8474 


35. Ni 31 -> A1 + A1 46 ->V 

-0.42323 

2.8518 

48.57 

36- Ni 30j32 -> V 

-0.27055 

2.8338 


37. Ni 30;31 -> V 

-0.35110 

2.8349 


38. Ni 30 -> A1 + Ni 31 -> V + Al 42 -> V 

-0.35407 

2.8383 


39. Ni 3 Q ; 3i -> A1 + Al 18 42 -> V 

-0.34075 

2.8442 

47.83 

40. Ni 2 6,30,34 _> v 

-0.34334 

2.8276 


41. Ni3o 5 3i 5 32 -> V 

-0.26174 

2.8263 


42. Ni 30j 34 ;35 -> V 

-0.24887 

2.8281 

47.06 

43. Ni 30; 3 lj54j 55 -> V 

-0.18313 

2.8218 


44. Ni 26; 27,29,30 "> V 

-0.14426 

2.8205 


45. Al 65 _ 68 -> V + Ni 53 _55 -> A1 

-0.10201 

2.8430 

46.27 

46. Ni 5 2 9 5 30,3 1,32,53 _> V 

-0.09769 

2.8070 


47. Ni 5j7j 29,3 1,53,55 “> V 

-0.06919 

2.7960 
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